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We study numerically the geometric properties of reduced super- 
symmetric non-compact SU(N) Yang-Mills integrals in D = 4 
dimensions, for N = 2,3,.. .,8. We show that in the range of 
^ large eigenvalues of the matrices A 11 , the original Z)-dimensional 

rotational symmetry is spontaneously broken and the dominat- 
ing field configurations become one-dimensional, as anticipated 
by studies of the underlying surface theory. We also discuss pos- 
sible implications of our results for the IKKT model. 

Introduction 

X 

Reduced supersymmetric Yang-Mills integrals appear in many areas of fun- 
damental physics ranging from QCD to string theory. An important issue 
is whether such reduced integrals are able to capture universal properties of 
the corresponding supersymmetric Yang-Mills theory in the large N limit, 
a question similar to the one posed by Eguchi and Kawai in the context of 
standard QCD |T|, 0]. Reduced supersymmetric Yang Mills integrals also 
enter calculations of the Witten index || of a wide class of supersymmetric 
quantum mechanical systems in D-brane physics 0, [J. Finally, in recent 
years a non-perturbative definition of string theory was proposed in terms of 
a ten-dimensional reduced supersymmetric Yang-Mills integral - the IKKT 
model - which is believed to reproduce, in an appropriate large N limit, 
amplitudes including all topological contributions of IIB critical strings H . 

As a supersymmetric model, the reduced Yang-Mills integral can be de- 
fined in D = 3,4,6, or 10 dimensions. As it turns out, many important 



1 



properties of the model, such as the singular behaviour of the field correla- 
tors, change with D in a systematic way, so that results gained for one value 
of D can be extrapolated to the other cases. We will in this paper focus 
on the four-dimensional model. In this case, the effective action, obtained 
by integrating out the fermionic degrees of freedom, is positive semi-definite 
which makes it well-suited to Monte Carlo studies. 

Analytically, the model - in any dimension - has so far been solved only 
for N = 2 |3], [5|, 0]. For larger N, there exist some results from approximate 



analytic methods M and numerical simulations || [TCJ, [11], [1^, |T3| . This made 



it possible to address questions related to the properties of the eigenvalue 



spectrum P, 10] as well as the scaling of physical quantities like the gyration 



radius, correlation functions, and Wilson loops 0, W\ [T3 



In this paper, we will study the geometrical properties of the model in 
more detail. Our motivation comes from earlier studies of the geometrical 
structure exhibited by the corresponding surface theory fI3| , |15f . This model 
was examined numerically in a discretized form, using the dynamical trian- 
gulation approach [JH], |T7| . The distribution of the gyration radius was found 
to have a power-like tail. 

For the matrix model one expects the same kind of large R behaviour. 
More precisely, the analytic solution for iV = 2 exhibits a tail p(R) ~ R~ 2D+5 , 
and it was conjectured from results of the one- loop approximation that this 
formula should also hold for any N > 2. Intriguingly, in the surface theory 
for D = 4 one observes p(R) ~ -R -3 , which agrees with the result from the 
matrix model. A mechanism was proposed to explain the appearance of this 
tail, relying on the emergence of essentially one-dimensional configurations 
called 'needles' or 'tubes'. Numerical simulations for D = 4 were able to 
show that such configurations do indeed dominate the large R part of the 
p(R) distribution. If one assumes the dominance of the 'tubes' for any D, 
one finds exactly the same formula as in the matrix model, p(R) ~ R~ 2D+5 . 

An obvious question that arises at this point is whether the tail that 
appears in the matrix model could be explained by a similar mechanism, 
i. e. whether the geometrical structure of large eigenvalue configurations 
also becomes one-dimensional. Trying to answer this question is the central 
aim of this paper. We will again use numerical simulations to examine this 
point. 

We begin the paper by shortly reviewing the matrix model of reduced 
Yang-Mills integrals. We then describe the observables that allow us to 
pick up the geometrical structure of a given configuration, and present the 
numerical results. We finish with a short discussion. 
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The model 



The zero-dimensional supersymmetric Yang-Mills integral is denned in the 
Euclidean sector by the partition function || 

Z = JdAdVdV e" 5 ^'* 1 , (1) 
where the action is given by 

S[A, *] = S B + S F = ~ Tr [A», A v } 2 + X - Tr ^T^A", *] . (2) 

The model is believed to provide a non-perturbative definition of string the- 
ory The N x N traceless Hermitean matrices A^, fi = 1, . . . , D, are a sort of 
quantum operators for the bosonic coordinates X M (r, a) — > A^ - of the string 
world-sheet in the D-dimensional target space. The fields ^ are fermionic 
matrices that transform as Majorana-Weyl spinors in D = 10 dimensions, 
and as Weyl spinors in D = 4. The IKKT model corresponds to the case 
D = 10. 

As in any quantum theory, one is interested in measuring correlation 
functions, the simplest of which are one-matrix correlators of the type 
{-^(TrAffi). Because the model is rotationally invariant, such correlators 
do not depend on the choice of A 11 . More specifically, for any matrix of the 
form A = n^A 11 , where is a unit vector, the correlators give identical 
results. 

One-matrix correlators can be expressed as moments of the distribution 
of eigenvalues of the A 11 , 

i(TM 2 )^= [ d\p(X)X 2k } (3) 

where 



For N = 2, the eigenvalue distributions can be found analytically for D = 4, 6 
and 10 dimensions 0, and were shown to have the large A behaviour 

p(X) ~ X-^ . (5) 

For higher iV the integrals have not been solved analytically. However, it 
was conjectured M, based on results of the one-loop approximation M, that 
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the large A part of the spectrum should be controlled by the same powers (||) 
independently of N. 

The power law (|5]) describes exactly the same behaviour as was found for 
the surface theory, in which case it was shown that the dominating configu- 



rations are one-dimensional |LJ] with scaleless fluctuations of the extension 
in the elongated direction. 

One way of extracting the dimensionality of the surface model is to mea- 
sure the correlation matrix 



fy/9 



c, u = j v ; r . (6) 



The trace of this matrix gives the square of the gyration radius, 

e>2 / v 7 ^ Y^ x l (7 , 

R -^ c »»- JVf ' ( ) 

We can simplify things by choosing a basis in which becomes diagonal, 

C tJ = ^S tJ rl (8) 

In this case, R 2 = r f- Also, the individual eigenvalues rf now describe the 
square extent of the system in the directions given by this particular basis. 

For large R, the distribution of R was found to behave as p(R) ~ -R -3 , 
just as in (|5]). Also, one of the eigenvalues rf was shown to become much 
larger than the others, reflecting the one- dimensional nature of the large R 
configurations in the surface theory [|17| . 



In analogy to this, we can now try to determine the dimensionality of the 
matrix model. We can define the quantum correlation matrix as 

= —TrA^At, (9) 

and the square of the gyration radius R 2 as 



R 2 



Given that we can freely rotate the N 1 as discussed above, we can again ar- 
range things such that we find a diagonal matrix Cy, and the gyration radius 
once more becomes R 2 = Y2i r i- The eigenvalues rf are all real numbers and 
can again be interpreted as the square extent of the system in the directions 
i, except that this time they also include quantum fluctuations. 
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The effective dimensionality can now be determined from the distribution 
of the rf. If d out of D eigenvalues can be shown to be much larger than the 
remaining ones, the system can be said to be effectively <i-dimensional. 

Let us order the eigenvalues by size, rf > rf +1 . We will study the be- 
haviour of rf as a function of R 2 : 

_ J 5(R 2 (A) - R 2 ) rf(A) 

2 r e -W>*l 5{R\A) -i? 2 ) 1 ' 

/,From what we found in the surface model, we expect for large R to see 
( r i)-R 2 ~ R 2 whereas the other (rf) R 2 should, in comparison, become negli- 
gible. It is convenient to introduce a new observable 

rl + rl + rl _ r\ 
T> R 2 _i R 2 [iZ) 

as a quantitative measure of this asymmetry. It tells us which fraction of 
the square extent of the system comes from the transverse directions. If 
the system does indeed become one-dimensional for large R, we should find 



Computational method 

Using the chiral representation of the gamma matrices and replacing the 
bispinors ^ by two-component spinors ip 




we can re-write the fermionic part of the action as 

S F = l -^ Maihbkn ^ kn (14) 
where the matrix M. is given by 

M mj , bkn = A^ kl f5 ni -A^f5 k3 (15) 

and the gamma matrices are 7^ = (— ia ,a). Here, space-time indices are 
denoted by Greek letters, matrix indices by i,j,k,..., and spinor indices 
by a,b,c,.... Each combination of sub-indices a,i,j forms a single index 
A = {aij} of the matrix M.ab- Since the possible values for the sub-indices 
are a = 1, 2 and i,j — l,...N,M is, for now, a 2N 2 x 2N 2 matrix. 



5 



However, as it stands M. has a pair of zero eigenvalues coming from a zero 
mode of the fermionic action, which is invariant under a change ip —>■ ip + el. 
This zero mode can be removed by omitting, in the partition function, the 
integration over one of the diagonal elements of the ip matrices, for example 
ipNN- This amounts to calculating the fermionic integral for a sub-matrix of 
M. in which the two rows and columns corresponding to % = j = N have 
been crossed out. Call the determinant of this sub- matrix detp(A). It is 
easily checked that detp(A) is always real and positive semi-definite. 

After integration over the fermions, the partition function becomes 

Z = J VA e" 5s(A) det F (A) (16) 
where the integration measure is given by 

va =u (ri rfRe K') rfim K)) (n^) j ( Tr ^) 

fJ. \i>j / \ i / 

The most practical way to deal with the global constraints TrA M = in 
(|T7j) is to just ignore them in the simulations, taking to be arbitrary 
Hermitean matrices. This creates a new zero mode in the bosonic sector, 
since the action is now invariant under the change A^ — > A^ + a M l, which 
manifests as a random walk of TyA^ in the simulations. However, this can 
easily be corrected by simply subtracting the trace from the matrices before 
each measuring step, 

K - A u - ^Tr^ . (18) 

Physically, this means that we always take measurements in a reference frame 
that is fixed to the system's center of mass. In the surface theory, the shift 
corresponds to replacing — > — X^ M . 

To update the fields we use a standard Metropolis algorithm. Each up- 
dating step consists of randomly choosing one element of each matrix N 1 and 
proposing a change by a random number taken from a uniform distribution 
in a range [— e, e], where e is adjusted so as to produce reasonable acceptance 
rates. We accept or reject the change according to the Metropolis criterion. 



Results 

A first, easy test of our program consists of measuring the average bosonic 
action, which we know from general scaling arguments should behave as 
(S B ) = 3/2(N 2 -1). 
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Figure 1: The distribution of eigenvalues from the theoretical prediction 
(solid line) and the numerical data (crosses). The scale is logarithmic on 
both axes. 

Secondly, we compared the distribution of eigenvalues of A^ 1 to the theo- 
retical formula for N = 2 



where U is the Kummer U-function. The numerical results are shown, to- 
gether with the theoretical curve, in figure |l|. 

For large A (large R 2 ) the quality of the numerical data is limited by two 
factors. For one thing, the number of data points in this range is very small as 
the tail, although long, contains only a small fraction of the whole distribution 
(^ 2.64% for |A| > 4). Thus, it takes a long time for the algorithm to produce 
reasonable statistics in this region. Secondly, the power of the tail is, in 
absolute values, not very large, which means that once the algorithm does 
enter the tail, it embarks on a random walk with long excursions that increase 
the autocorrelation time enormously. Despite these limitations, however, the 
figure shows quite good agreement between the theoretical curve and the 
simulation data even within the tail. 

Generally, it is known that for the reasons just given it is extremely 
difficult to deal with power-like fall-offs in numerical simulations. To improve 




(19) 



the quality of the data in this region, we can set a lower limit R^in on 
to prevent the algorithm from going to the bulk of the distribution, where it 
would spend almost all of its time otherwise. The price to pay for such a cut- 
off is a decrease in the acceptance rate of the algorithm as it frequently tries 
to push through this boundary and go back to the main part. In practice, 
however, the drop in the acceptance rate turns out to be not too severe. For 
example, with N = 4 and e = 0.1 the acceptance rate decreases from 74% to 
24% after imposing a lower barrier of R 2 min = 120. 

Similarly, to prevent the algorithm from making too long excursions into 
the comparatively flat tail of the distribution we can introduce an upper limit 
R 2 max - As expected, this drastically reduces the autocorrelation time. 

We measured the asymmetry parameter rj and the distribution of the 
squared gyration radius R? for N — 2, . . . , 8. For each value of N we per- 
formed 10 5 measurements (or 4 x 10 4 in the case of N — 8), each separated by 
100 sweeps, where a sweep encompasses N 2 Metropolis trials. As an example, 
for iV = 4, we found an integrated autocorrelation time for R 2 , calculated in 
units of 100 sweeps, of t R 2 = 130(10). 

Let us start with a discussion of the results from the N = 2 case. From 
the analytic solution (P, W), we find for the distribution of the eigenvalues 



of 



p {r) „ ( rir2r3 )«( r 2_ r 2)( r 2_ r 2^ r 2_ r 2^-2(rfri+r|ri+r|r?) _ ^ 



All eigenvalues r^, % > 3, are identically zero, i. e. even for small R 2 the sys- 
tem is only three-dimensional independently of D. The exponent a depends 
on D as a = 2D — 5. The same formula, with a = D — 3, also describes the 
distribution of eigenvalues of the purely bosonic model. 

We can use equation ( |2"0"D to explain the large R behaviour of the model, 
where in this case R 2 = r\ + r\ + r| (because r\ = 0). When moving to 
large values of R 2 , all configurations will be exponentially suppressed, except 
for those that keep the exponent r\r\ + r\r\ + r\r\ approximately constant. 
Given that R 2 can only grow large if at least one of the eigenvalues does so 
as well, this can be achieved only if r\ ~ R 2 , r%3 ~ -jp- As a consequence we 
can, in (|20|) , approximately replace r 1 — > R, r 2 ^ — > R" 1 and f dr 2t 3 — > R~ x . 
This gives us the distribution of R in the large R limit, 

p(R) ~ R~ a . (21) 
Since we have -|# ~ ^4 (see figure 0), we find for the asymmetry parameter 

"~Wr- (22) 
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Figure 2: The eigenvalues r 2 , rf, and rf, normalized with R 2 , as functions 
of R 2 for N = 2. Note that the data has been smoothed so that each point 
in the figure actually represents a pair (R 2 ,f]), with the average taken over 
n = 100 successive values of R 2 . 

i. e. 7] goes to zero for large R 2 and the quantum system becomes one- 
dimensional. 

To summarize for N = 2: The flat directions of the bosonic part of the 
action, which correspond to constant values of the exponent in (Iffi), lead 
to a large R behaviour of the eigenvalues r\ ~ R 2 , r| 3 ~ -j^. This is true 
independently of the presence or absence of fermions in the theory. However, 
the addition of fermions does strongly influence the power of the probability 
distribution p(R). 

Note that the flat directions being one-dimensional is a consequence of 
the particular form of the bosonic action; we can create systems of higher 
dimensions in the large R range by choosing a different action. For example, 
the action Tr [A M , A u ] [A u , A a ] [A a , A M ] has the same symmetries as TrfA^, A u } 2 , 
but some of its flat directions are now two-dimensional. This can be seen 
from the eigenvalue distribution, which in this case reads 

p(r) ~ (r 1 r 2 r3) a (r?-r 2 2 )(r 2 2 -r3 2 )(r 3 2 -r2)e- 24 ^^'-3 . (23) 

For large R, one can keep the exponent r\r\r\ constant by choosing r\ ~ R, 
r 2 ~ R, and r 3 ~ 1/R, which obviously describes a two-dimensional disc. In 
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Figure 3: The plot shows (upper curve) and r\ = T2 ^ 3 2 Ts (lower curve) 
as functions of R 2 , for N = 4. The data has been smoothed as explained 
for figure |2|. The plot combines two data sets, one from an unconstrained 
simulation (the points for R 2 < 100) and one from a simulation with a lower 
cut-off R 2 min = 120. The solid line is the best power-like fit to the data from 
the unconstrained simulation. 

this case rj(R) goes toward a constant for large R. This type of asymmetry 
could be traced by a higher dimensional counter-term of rj, rj 2 = 1 — Tl ^ 2 , 
which would vanish for large R. 

Coming back to the standard model, we can now use Monte Carlo simu- 
lations to repeat the same analysis for N > 3. Figure [3] shows the numerical 
data for rj measured as a function of R 2 for N = 4. As for N = 2, for large R 2 , 
7] vanishes - in other words, the main contribution to the square extent comes 
from the elongated direction of what is thus an essentially one-dimensional 
system. 

We fitted the data to a power-law function rj = a(R 2 ) b , and found for 
iV = 4 as the best fit a = 1.378(7) and b = -0.690(4). To test the stability 
of this result in the region of large R 2 , we repeated the simulation with the 
inclusion of a lower bound R 2 min = 120 as discussed in the previous section, 
to see whether the validity of the fit goes beyond the range of R 2 generated 
in the first run. The results are collected in figure and do indeed show 
good agreement with the results of the first simulation. 
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Finally, we collect the results of the best power-law fits to 77 for iV = 
3, 4, 5, 6, 8 in the following table. For comparison, we also include the results 
for the surface model with = 8, 12, 20, 28 triangles. 



N 


bmatrix 


tit 


bsurface 


3 


-0.825(2) 


8 


-1.604(16) 


4 


-0.690(4) 


12 


-1.463(19) 


5 


-0.606(4) 


20 


-1.412(15) 


6 


-0.562(5) 


28 


-1.383(20) 


8 


-0.503(29) 







In both the matrix and the surface model, the exponent b changes with N. If 
we assume the finite size corrections to be of order that is b(N) = b oc + -^, 
we can try to estimate from a fit to the measured values. For the matrix 
model, this leads to b^ = —0.35(1), c = —1.37(4). We also see that b^ 
stays negative and far from zero if we assume a different kind of finite size 
correction, like or Thus, the data suggests that the dominance of 
one-dimensional configurations in the large R part of the spectrum should 
also be present in the limit of N — > 00. To answer this question conclusively, 
however, one should extend the simulations to larger values of N. 

Another related issue is the scaling of the gyration radius with N. From 
results of the one- loop approximation [§], one expects the maximum of the 
gyration radius distribution to scale as ■sj R max (N) ~ N 1 ^, corresponding 
to the Hausdorff dimension du — 4 of branched polymers which appear 
naturally in this approximation. Indeed, if we fit our numerical results for 
the surface model to the formula \J R max {N) = aN b , we find a = 0.573(4) 
and b = 0.234(4), which is close to 1/4 even though the surfaces we studied 
are quite small = 12, ... , 60). Alternatively, we can do the same using 
another quantity R, which was defined in |iT| as an alternative estimate of 
the system extent : 

i 

Here, the best power-like fit gives a = 0.579(6) and b = 0.204(4), which 
does not differ too much from the values obtained for the generic definition 
of the gyration radius. Presumably, for larger values of A^ we should see b 
converging to 1/4 for either definition. 

In the matrix model, on the other hand, the situation looks quite different 
in the simulated range of N. Namely, in this case the best fit to ^R 2 max {N) = 
aN b is a = 1.03(4) and b = 0.81(6), which means d H = \ = 1.23(9). This is 
far away from the result of the one-loop approximation, d H = 4. However, 
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one should note that the simulated systems are rather small for branched 



polymers to fully develop. Indeed, in [|Tl|] a semi-classical analogon of R 
constructed from commuting coordinates was studied for values of N in the 
range N — 16, . . . ,48, and its average was estimated to behave as (R sc )n ~ 
iV 1//4 , which should also mean (VR 2 ~) N ~ iV 1//4 . 

Let us close this section with a remark about the definition of dimen- 
sionality of the system. There are two natural possibilities : the Hausdorff 
dimension d^, and the effective dimension coming from the principal compo- 
nent analysis of the correlation matrix C^ u . The former describes the relation 
between the system size and its average extent, whereas the latter gives the 
minimal dimensionality of the subspace to which we can restrict our descrip- 
tion of the system without neglecting important degrees of freedom. The 
two may differ in general, and it is an important physical question which 
definition should be used for any given problem 



Conclusions 

We have shown that the underlying geometry that can be associated with 
the field configurations of reduced supersymmetric Yang-Mills integrals in 
D = 4 dimensions shrinks in the limit of large R 2 = Tr A 2 to an essentially 
one-dimensional tube. The original rotational symmetry of the action is 
spontaneously broken to the direction of the tube. The same mechanism has 
already been observed in numerical simulations of the corresponding surface 
model. 

The origin of the power-like behaviour boils down to the existence of flat 
directions in the bosonic part of the action. The tubes correspond to field 
configurations that expand along the valleys of these flat directions, where 
they do not have to pay the usual 'exponential price', but rather one that 
is only power-like. The exact value of this power comes from the prefactor, 
which depends on the dimensionality of the problem. However, it seems 
that the power does not change with N. For D = 4, the distribution of the 
gyration radius was determined as p(R) ~ R~ 3 . 

The analysis presented in this paper should be applicable to the six- 
and ten-dimensional cases as well, where one expects a distribution p(R) ~ 
R~ 2D+5 . In the large N limit we expect the estimator of the linear system 
extent as defined by Rk = (R 2k )^ to pick up the singularity for k large 
enough, which may lead to a different scaling with N and thus a different 
Hausdorff dimension. A similar thing happens for example with Levy random 
walks, which are known to have a Hausdorff dimension smaller than two 



19. For the IKKT model in particular, an eigenvalue distribution with a 
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power-like tail as the one suggested here would mean that we should expect 
higher-order correlators starting with k > 7 to effectively sample the tail of 
the distribution, which is dominated by one-dimensional configurations. 

Acknowledgments 

We thank Jan Ambj0rn, Piotr Bialas, Jun Nishimura and Matthias Stau- 
dacher for discussions. The work was partially supported by the EC IHP 
network HPRN-CT-1999-000161 and the KBN grant 2P03B00814. 

References 

T. Eguchi and H. Kawai, Phys. Rev. Lett. 48 (1982) 1063. 
[2] D. J. Gross and Y. Kitazawa, Nucl. Phys. B206 (1982) 440. 
[3] A. V. Smilga, Nucl. Phys. B266 (1986) 45. 
[4] P. Yi, Nucl. Phys. B505 (1997) 307. 

[5] S. Sethi and M. Stern, Comm. Math. Phys. 194 (1998) 675. 



[6 
[7 



[9 

[10 
[11 

[12 
[13 



N. Ishibashi, H. Kawai, Y. Kitazawa and A. Tsuchiya, Nucl. Phys. 
B498 (1997) 467. 

W. Krauth, J. Plefka and M. Staudacher, Class. Quant. Grav. 17 (2000) 
1171-1179. 

H. Aoki, S. Iso, H. Kawai, Y. Kitazawa and A. Tada, Prog. Theor. 
Phys. 99 (1998) 713. 

W. Krauth, H. Nicolai and M. Staudacher, Phys. Lett. B431 (1998) 
31. 

W. Krauth and M. Staudacher, Phys. Lett. B453 (1999) 253. 

J. Ambj0rn, K. N. Anagnostopoulos, W. Bietenholz, T. Hotta and J. 
Nishimura, JHEP 0007 (2000) 013. 

J. Ambj0rn, K. N. Anagnostopoulos, W. Bietenholz, T. Hotta and J. 
Nishimura, |hep-lat/000903C| . 



J. Ambj0rn, K. N. Anagnostopoulos, W. Bietenholz, T. Hotta and J. 
Nishimura, JHEP 0007 (2000) 011. 



13 



[14] I. Bars, Phys. Lett. B245 (1990) 35. 

[15] M. Fukuma, H. Kawai, Y. Kitazawa and A. Tsuchiya, Nucl. Phys. B510 
(1998) 158. 

[16] S. Oda and T. Yukawa, Prog. Theor. Phys. 102 (1999) 215. 

[17] P. Bialas, Z. Burda, B. Petersson and J. Tabaczek, Nucl. Phys. B591 
(2000) 1. 

[18] J. Nishimura and G. Vernizzi, JHEP 0004 (2000) 015. 

[19] J.-P. Bouchaud and A. Georges, Phys. Rep. 195 (1990) 127. 



14 



